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We study the mechanism of decay of a topological (winding-number) excitation due to finite- 
size effects in a two-dimensional valence-bond solid state, realized in an S = 1/2 spin model ( J-Q 
model) with six-spin interactions and studied using projector Monte Carlo simulations in the valence 
bond basis. A topological excitation with winding number |Wj > 0 contains domain walls, which 
are unstable due to the emergence of long valence bonds in the wave function, unlike in effective 
descriptions with the quantum dimer model (which by construction includes only short bonds). We 
find that the life time of the winding number in imaginary time, which is directly accessible in the 
simulations, diverges as a power of the system length L. The energy can be computed within this 
time (i.e., it converges toward a “quasi-eigenvalue” before the winding number decays) and agrees 
for large L with the domain-wall energy computed in an open lattice with boundary modifications 
enforcing a domain wall. Constructing a simplified two-state model which can be solved in real 
and imaginary time, and using the imaginary-time behavior from the simulations as input, we find 
that the real-time decay rate out of the initial winding sector is exponentially small in L. Thus, 
the winding number rapidly becomes a well-defined conserved quantum number for large systems, 
supporting the conclusions reached by computing the energy quasi-eigenvalues. Including Heisenberg 
exchange interactions which brings the system to a quantum-critical point separating the valence- 
bond solid from an antiferromagnetic ground state (the putative “deconfined” quantum-critical 
point), we can also converge the domain wall energy here and find that it decays as a power-law 
of the system size. Thus, the winding number is an emergent quantum number also at the critical 
point, with all winding number sectors becoming degenerate in the thermodynamic limit. This 
supports the description of the critical point in terms of a U(l) gauge-field theory. 

PACS numbers: 75.10.Kt, 75.10.Jm, 75.40.Mg, 75.40.-s 


I. INTRODUCTION 

Systems with topological order are characterized by 
unconventional quantum numbers labeling degenerate 
ground states, the number of which depends on bound¬ 
ary conditions. Except for models where conservation of 
topological numbers is ensured by construction, such as 
Kitaev’s toric code^ in general a topological quantum 
number is only emergent in the limit of infinite system 
size. For a finite system the life time of a state prepared 
with a fixed topological number is finite, and, thus, the 
levels within the ground-state manifold are split. One 
would normally expect the life time to be exponentially 
long in the system size, as in the case of an ordered state 
with a broken discrete symmetry (e.g., an Ising ferro- 
magnet in a weak transverse field). In practice, e.g., 
when considering the design of topologically protected 
qubits*^— the life time due to finite system size may then 
not play an important role if the qubit is sufficiently large. 
It is still useful to investigate quantitatively the mecha¬ 
nism of these instabilities due to finite size in various 
topological states. The emergence of topological conser¬ 
vation laws is also of interest in descriptions of quantum 
matter using simplified Hamiltonians and quantum field 
theories where they are conserved by construction, e.g., 
in dimer models and field theories derived from them. 7 ' 8 

Topological quantum numbers are often discussed in 


the context of quantum spin liquids* — but can also ap¬ 
pear in systems with long-range order. Here we inves¬ 
tigate a quantum spin 5 = 1/2 Hamiltonian which has 
a valence bond (VB) ordered (spontaneously dimerized) 
ground state and whose excitations containing domain 
walls can be classified by a winding number (which es¬ 
sentially counts the number of domain walls in the sys¬ 
tem). In this case the ground state of a large system 
is in one winding number sector and the other sectors 
are at higher energy. Our interest here is to study quan¬ 
titatively the instabilities of the domain walls and the 
topological mechanisms (changes of the winding num¬ 
ber) responsible for their decay in finite systems. This 
issue is of particular interest when the VB solid order is 
weakened by the introduction of interactions that eventu¬ 
ally completely destroy the order and bring about other 
phases. In this process the winding numbers should also 
become unstable in the thermodynamic limit. In our 
studies discussed here, we use the so-called J-Q model 
on the square lattice ; 13 ’ 14 where J is the standard anti¬ 
ferromagnetic Heisenberg exchange and Q a multi-spin 
(here six-spin) interactions favoring VB order. 

The singlet sector of the SU(2) invariant spin system 
considered here is amenable to a description in the VB 
basisji^— where a state of an even number of spins N is 
a superposition of product states containing N/2 singlet 
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pairs (VBs), for which we use the notation 

(«, b) = (UU - UW/V2. (1) 

While the winding number is strictly conserved within 
a restricted basis of short VBs (bond lengths < L/4, 
where L is the system length, as we will explain in 
detail below) ^ with all bond lengths included, as re¬ 
quired for the basis to be complete (whence the basis in 
fact becomes over-complete), the winding number is no 
longer conserved, and, thus, domain walls can decay due 
to topological fluctuations. We here study such decay 
within the imaginary-time dynamics accessible in pro¬ 
jector quantum Monte Carlo (PQMC) simulations in the 
VB basis ; 20 ! 21 

We find that an initial state with domain walls decays 
to the ground state with no domain walls through tran¬ 
sitional states with a high density of long VBs. As a re¬ 
sult, the life time is only growing with the system size as a 
power law. The analytic continuation between imaginary 
and real time is very complex for the non-equilibrium 
situation considered here, and we cannot translate the 
imaginary-time behavior rigorously to real-time evolu¬ 
tion, where an isolated system should thernralize at con¬ 
stant energy instead of decaying to the ground state. To 
gain some insights into the relevant real-time scale cor¬ 
responding to the power-law divergent imaginary-time 
scale of the winding number, we consider an effective 
two-state model in combination with the PQMC data. 
Based on this approximation an exponentially long life 
time in real time appears plausible, both in the VB solid 
phase and at the critical point separating it from an an- 
tiferromagnetic ground state. 

In addition to the life-time of the winding number, we 
also discuss the quasi-eigen energies of the winding states 
in the VB solid and at the critical point. Comparisons 
with domain-wall calculations in boundary-modified sys¬ 
tems with forced domains confirm the domain-wall na¬ 
ture of the states with non-zero winding number. 

In Sec. mi we will discuss the winding number in detail 
and explain our methods to investigate its fluctuations 
in PQMC simulations. Results of such simulations in the 
VB solid phase of the J-Q model are presented in Sec. m 
In Sec. m we introduce the effective two-state model 
and study the relationship between imaginary-time de¬ 
cay to the ground state and transition rates in real time. 
We consider scaling in a quantum-critical system and the 
eventual instability of the winding number sectors in the 
antiferromagnetic phase in Sec.fV] In sec. IVIl we summa¬ 
rize and further discuss our results and their implications. 

II. DIMERS, VALENCE BONDS, AND 

WINDING NUMBERS 

To more precisely introduce the concepts and mech¬ 
anisms to be discussed below, consider first the wind¬ 
ing number W = (w x ,w v ) of a classical close-packed 
dimer model on the square lattice. A dimer connects two 


nearest-neighbor sites, one on sublattice A and one on 
B 1 with the A and B sites forming a checker-board pat¬ 
tern. W can be defined by assigning a direction (arrow) 
A—>B for each dimer. Superimposing any such configu¬ 
ration onto a reference configuration with B—»A dimers 
forming (by convention) horizontal columns, closed loops 
form and w x , w y correspond to the x and y currents nor¬ 
malized by the system length L. 

The classical dimer configurations are the basis states 
of the quantum dimer model (QDM)*- 22 ~— 4 which in the 
simplest case has a diagonal term counting the number 
of pairs of parallel bonds and an off-diagonal term which 
can rotate such “flippable pairs” by 90°. The off-diagonal 
terms being local, they cannot change W, which, thus, is 
a good quantum number of the QDM on a periodic lat¬ 
tice. The Hilbert space consists of ~ L 2 winding number 
sectors. On a non-bipartite lattice, e.g., the triangular 
lattice, or in an extended square-lattice model includ¬ 
ing also bonds connecting next-nearest neighbors^ some 
winding numbers mix and the sectors are reduced down 
to even and odd ones / 9 ' 26 labeled, e.g., by w x ,w y = ±1. 
On the torus there are thus four sectors of conserved W. 

Now consider S = 1/2 spins. Any total spin sin¬ 
glet can be written as a superposition of tilings with 
VBs, and if the system is bipartite one can restrict the 
bonds to only connect sites on different sublattices. If 
the subscripts a and b in Eq. © correspond to the sub¬ 
lattices as a £ A, b £ B, the sign of each component of 
|u) = |(ai, b\) ■ ■ ■ (a^v/ 2 ; £>jv/ 2 )) i n the t,4- basis conforms 
with Marshall’s sign rule for the ground state of a bipar¬ 
tite system. Such a state can be expanded in VB states, 

1*0 > = £/», (2) 

V 

with positive-definite expansion coefficients f v (where 
v labels the bond configurations, v £ {1,..., (V/2)!}, 
which map into permutations of the N/2 A sites con¬ 
nected to N/2 B sites/ 16 ’ 17 The VB basis states are 
non-orthogonal, with (u 2 |vi) = 2 ni2-Ar / 2 , m 2 being 
the number of loops in the transition graph of the VB 
configurations ^ 6 - 1 ' 

The insight that many quantum states of spins can, 
to a good approximation, be expressed with short VBs 
motivated the introduction of the QDM as a class of ef¬ 
fective models to describe some quantum spin systems, 
spin liquids in particular ^2 r 24 i 27 ~ — To judge whether 
the QDM indeed provides a good description in a given 
case, one has to consider the role of long bonds, the 
non-orthogonality of the VB basis, and the interactions 
included in the QDM (for which there is a systematic 
schenr o 31 ! 32 ). In general it is not possible to rigorously 
prove the validity of the QDM description, other than by 
careful comparisons of numerical results when available. 

The QDM should provide an excellent approximation 
for a VB solid (crystal), where predominantly short VBs 
form a regular density pattern. Such ground states 
are hosted by many quantum spin models with frus¬ 
trated exchange interactions in parts of their parameter 
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FIG. 1. (Color online) VB configurations with winding num¬ 
ber (0,0) in (a) and (1,0) in (b) and (c), with (c) obtained 
from (b) by locally rotating pairs of dimers. Configuration (a) 
serves as the reference configuration for defining the winding 
numbers. These VB states are all used as initial states in the 
PQMC simulations discussed in this paper. 


spacesj22r— and also by models with certain multi-spin 
couplings beyond the pair exchange d 3 ' 14 ! 38 ' 39 One of the 
latter is the J-Q model^ where J is the antiferromag¬ 
netic Heisenberg exchange, which we can write as 

2 N 

H J = - J ^2 > ( 3 ) 

6=1 

where Pij is a singlet projector, 

Pij ' = (1/4-Si- S,-), (4) 

and b labels the links connecting nearest neighbors 
i(b),j(b) on the square lattice with N = L 2 sites. For 
the Q-term we use the six-spin variant^ 

2N 

Hq = Pi( c )j(c)Pk(c)l(c)Pm(c)n(c), (b) 

c=1 

where the cells c contain three nearest-neighbor site pairs 
(■ ij){kl){mn ) forming 2x3 and 3x2 plaquettes. The 
model H(J , Q) = Hj + Hq has a strongly ordered colum¬ 
nar VB solid ground state for large Q/J 42 We will here 
first study H = Hq as a stand-alone models and later 
tune Q/J to the critical point at which the VB solid order 
vanishes and antiferromagnetic order sets in. 

Performing PQMC simulations in the over-complete 
bipartite VB basis, we have access to the winding num¬ 
ber and our aim is to study its stability as a function 
of the system size. A columnar arrangement of short 
VBs has winding number W = (0,0) and any other W 
corresponds to the presence of domain walls, as illus¬ 
trated in Fig. [TJ Domain walls have been studied within 
QDMs^ including quantum phase transitions driven by 
increasing winding number (density of domain walls) i± 
Here our motivation is different, but it is useful to have 
a VB solid of the QDM as a reference point, where the 
winding numbers characterizing states with domain walls 
are fully conserved by construction. 

In this study we use the imaginary-time Schrodinger 
evolution operator 


U{t) = e~ rH , ( 6 ) 



FIG. 2. (Color online) Illustration of the winding number 
and its non-conserved property for an 8-site chain. The open 
and solid circles indicate sublattices A and B, respectively, 
and the lines and arcs with arrows are the VBs. In (a) the 
reference configuration for defining W is shown, and in (b) it 
has been superimposed with a state with W = 1. In (c) a 
singlet projector has acted on the W = 1 state, at sites 3 and 
4 from the left, and in (d) a subsequent action on sites 5 and 
6 has lead to a change of the winding number to W = 0 (the 
bonds lost in the processes are shown with dashes). 


and apply it in PQMC simulations^ to an initial state 
|^ / o(w ; x)) with only short (nearest-neighbor) bonds and 
winding number [w x , 0). A w x ^ 0 state is obtained 
from a perfect columnar state by shifting an odd number 
of rows, as illustrated in Fig. |Tj for w x = 1. We have 
U(t —> oo)|\l/o(w x )) —> |0) (neglecting an unimportant 
normalization), where |0) is the ground state, which is 
dominated by the w x = 0 sector; for L —> oo the prob¬ 
abilities approach P(w x = 0) = 1 and P{w x = 1) = 0 
exponentially with increasing L. 

In the program implementation of the above PQMC 
approach, we Taylor expand the operator © to all con¬ 
tributing orders, as in the stochastic series expansion 
(SSE) method.- 2,43 Thus, each configuration consists of 
a series of n oc Nt operators which are sampled out of 
the two- and six-spin terms in Eqs. m and (151) . In the 
model with J = 0 and Q = 1 that we will study first, 
there is, thus, a total of 3n singlet projectors acting on 
the initial state for a term of order n, while with both 
J > 0 and Q > 0 there are both two- and six-spin oper¬ 
ators present and the total number of singlet projectors 
is between n and 3n. 

To see how the presence of long bonds leads to non¬ 
conservation of W we discuss a one-dinrensional example 
with 8 spins, illustrated in Fig. [2] The reference config¬ 
uration used to define W £ {0,1}, which by definition 
itself has W = 0, is shown in (a), while (b) shows its 
transition graph with the W = 1 short-bond state. In 
(c) a singlet projector has acted on the W = 1 state and 
reconfigured two bonds, leading to a new bond of length 
1 and one of length 3. The winding number remains 
at W = 1. A second operation in (d) leads to a bond 
which in a longer chain would have length 5, but in this 
L — 8 periodic system has length L — 5 = 3 by definition 
of the VB basis. Then the winding number changes to 
W = 0. This is also in accord with an examination of the 
bond pattern, which now has the short bonds shifted by 
one lattice spacing relative to those in the original bond 
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configuration. Quite generally, if the two bonds affected 
by a singlet projection span more than half the system 
length, then the winding number changes in the process, 
i.e., the minimum bond length required for W to not be 
conserved is L /4 for L being a multiple of 4. We study 
also L given by an odd multiple of 2, for which there are 
only small, trivial differences from the above in how the 
winding numbers change. 


III. RESULTS IN THE STRONGLY-ORDERED 
VALENCE-BOND SOLID STATE 



If the winding number is conserved for L —> oo, 
which we expect, then in this limit the projected state 
U (t) | 4/q (re3,)} should evolve toward the lowest eigenstate 
with winding number w x . For finite L, one can ex¬ 
pect there to be a “quasi-eigenstate” toward which the 
state evolves before w x typically decays to lower val¬ 
ues (with some fluctuations possible also to higher w x ) 
and the ground state is obtained. Sampling the norm 
(4to(w x )|/7(T)C/(r)|'ho(w a; )}, we can use the estimator 
(H) = — (n}/(2r), in analogy with the SSE method . — 
We classify a PQMC configuration as having a conserved 
winding number if propagation with U 2 from both the 
left and the right maintains w x in its initial sector af¬ 
ter each step. This way, we can compute the energy in 
different sectors, and also in the ground state for large 
enough r. We have reported preliminary results for the 
domain wall energies in different winding sectors based 
on a slightly different procedure . 44 ' 45 

In Fig. [3] we show results for the quantity 


k(w x ,L) 


(H) Wx - (H )o 
4 Lw x 


(7) 


which can be interpreted as a domain wall energy per 
unit length when converged. We previously computed 
the domain wall energy based on systems with edge mod¬ 
ifications favoring VB ordering in such a way that a single 
domain wall of the type in Fig. 0b) is present or absent.— 
Such a domain wall can be classified as having a twist an¬ 
gle 4> = 7T ) 40,46 while a domain wall between horizontal 
and vertical VB solids, as in in Fig. [He), has <)> = 7t/2 
(and the ir wall of course consists of two separate 7t/2 
walls). In periodic systems the total VB twist angle is 
27 tw x . For the open systems we therefore define n corre¬ 
sponding to Eq. 0 by dividing the energy difference of 
systems with and without domain walls by 2L<//7 t. As 
can be seen in Fig. 01 the results of different calculations 
give consistent results for k when L —> oo, confirming the 
above arguments. 

Simulations for large periodic systems suffer from er- 
godicity problems (long equilibration times) due to which 
the system may stay for a very long time in the initial w x 
sector, even at projection times r where w x = 0 should 
dominate (as judged by the behavior for smaller systems). 
In the energy calculation this is an advantage, as it al¬ 
lows us to converge very well to the w x ^ 0 lowest quasi¬ 


FIG. 3. (Color online) Domain wall energies per unit length, 
normalized to a single domain wall with twist angle p = tt/2. 
Results are shown for periodic systems in the sectors w x = 1,2 
and for systems with open-boundary modifications enforcing 
domain walls. For the periodic systems the projection time 
was t = L/8. The angle 4> indicated is the total twist of the 
VB order parameter when going around a periodic system or 
across an open system. 
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FIG. 4. (Color online) The probability of a state projected 
out of an initial state with w x = 1 to remain in that winding 
number sector at time r. 


eigenstates, in a way similar to measuring the energy of 
a meta-stable state in classical Monte Carlo simulations. 
Using t = L/8 and an initial state with w x = 1, the cal¬ 
culations for L up to 32 in Fig. 0showed fluctuating w x , 
while runs for larger L typically stayed locked at w x = 1. 

The imaginary-time life time of the state evolved from 
|iko(w x )) can be defined, e.g., as the r at which the 
probability of remaining in the initial w x sector is 1/2. 
However, because of the aforementioned ergodicity prob¬ 
lems this definition is practically useful only for relatively 
small system sizes. We have therefore developed an al¬ 
ternative approach, by using a long total projection time 
/3 = MA t , writing U(/3) as U M ( A r ), and individually 
Taylor-expanding each of these factors in the PQMC sim¬ 
ulations. Then we can monitor the winding number in 
the state propagated from the left and from the right 
with the full operator string (i.e., the concatenation of 
the M individual strings) and in both cases there will be 
some time “slice” r = iA T (with i counted from the bra 
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FIG. 5. (Color online) Life time, defined as the projection 
time r = T 1/2 at which P(w x ) = 1/2 (based on interpolation 
of data such as those in Fig. II , versus the system size on a log- 
log scale, (a) shows the behavior for w x = 1 with two different 
initial states; generalizations of those depicted in Fig. Eb) 
and (c). The lines show the forms ~ L°' j5 and ^ T 0 -' drawn 
through the large-L data obtained with initial states Figs. 1 
(b) and 1 (c), respectively. In (b) the life times for w x = 1,2, 
and L/ 2 are shown [using initial states of type Fig. Eb)], 
along with the form ~ L 0,55 drawn through the large-L data. 

or ket initial state) at which the winding number changes 
from, say w x = 1 to 0. Once a system has equilibrated 
and /? is large enough (so that the energy has converged), 
we can measure the probability of staying in the initial 
w x sector as a function of r. This method has much less 
severe autocorrelation problems once a decay to w x = 0 
somewhere in the time space has occurred. 

Some typical results for w x = 1 are shown in Fig. E 
and the resulting life time Ti/ 2 (L), defined as the r for 
which the probability P(w x = 1) of the state to stay in 
the initial sector equals 1/2, is graphed in Fig. [5] While 
this definition of t\/ 2 is not, strictly speaking, based on a 
bona fide quantum mechanical expectation value, it nev¬ 
ertheless gives a life time of the same order as the original 
definition proposed earlier (and this can also be expected 
based on formal considerations of the dynamics generated 
by operator products^). For the system sizes where we 
have data available from both approaches, the original 
definition gives a somewhat larger value but similar size 
dependence. As seen in Fig. [3(a), after a cross-over be¬ 
havior for small L , the life time grows asymptotically as 
a power law L a . 

While the power-law behavior appears to be robust 
to variations in the initial state, the value of the expo¬ 
nent may not be. Using two different initial states of 
the types (b) and (c) depicted in Fig. [3 two different 



FIG. 6 . (Color online) Bond-length probability for bonds of 
shape (a:, 0 ) at projection time r = tj / 2 for different system 
sizes when the starting state has winding number w x = 1 [of 
the type of the simple domain-wall state depicted in Fig.Eb)] • 
For reference, results are also shown for the VB solid ground 
state (W = 0 sector). 

exponents a are obtained. The initial states both have 
w x = 1 , but (c) explicitly implements ir/2 domain walls, 
while the state in (b) has sharp n domain walls. The ir 
domain walls will split up into n/2 domain walls in the 
course of imaginary-time projection, and by implement¬ 
ing them from the outset the initial state is closer to the 
eventual w x = 1 quasi-eigenstate. The range of system 
sizes for which the power-law can be approximately fitted 
in Fig. E a ) is rather small, and the exponents may still 
drift for larger systems. We can therefore not exclude 
that the exponents are asymptotically the same, though 
it also appears reasonable that state (b) is shorter-lived 
by a power of L, due to it being further away from the 
quasi-eigen state because of the wrong kind of initial do¬ 
main walls imposed. 

In Fig. [3(b) we compare results for w x = 1, w x = 2, 
and the extreme case of w x = L/2. In all cases the initial 
state was of the type in Fig. Eb), with only horizontal 
bonds shifted to achieve the different winding numbers. 
The life time for w x > 1 is also defined based on re¬ 
sults such as those in Fig. E as the probability of still 
remaining in the original winding sector being 1/2. The 
initial event of decay is almost always into a state with 
w x one unit smaller than the initial value. In Fig.Eb) we 
have fitted the results for all w x and large L to the same 
power law L 0 - 55 , though, again, there are considerable 
uncertainties in the exponents and we cannot conclude 
positively that they really are the same in all cases. The 
prefactor of the power-law decreases with increasing w x , 
but, interestingly, even for w x = L/2 it is not that much 
smaller than at w x = 1 . 

With only short VBs the life time as defined above 
would be infinite. With the probability of long bonds 
decaying exponentially with the bond length in both the 
initial and final states (which we have confirmed), one 
might naively have expected an exponentially long life 
time. The reason for the much shorter life-time is that 
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FIG. 7. (Color online) Log-log plot of the total probability of 
VBs of length > L /4 in the ground state and in the transient 
w x = 1 states at n/2 versus the system size, obtained under 
the same conditions as the results in Fig. [6] The (red) straight 
line has slope —2, while the curve (black) shows the form 
Pw oc exp(— aL b ) with b ~ 0.7. 


populated. In principle, our T 1/2 computed in the pre¬ 
vious section should in some way (through analytic con¬ 
tinuation) be related to this time scale. While in some 
cases it is understood how to relate real and imaginary- 
time evolution^ in the present case there is no known 
way to do this based on numerical imaginary-time data. 

To gain some insights into the real-time dynamics cor¬ 
responding to the power-law increase in the imaginary 
life time with L , we consider an approximate but illu¬ 
minating simple effective two-state model. We consider 
unperturbed states | 4-) and | f) corresponding to the 
w x = 0,1 sectors of the VB solid for large L , with ener¬ 
gies E± = —e and E^ = e. To mimic the decay of a state 
initially in the higher-energy w x = 1 sector to wii, = 0, 
we consider a perturbation by an off-diagonal matrix el¬ 
ement x <C e, i.e., the effective Hamiltonian 

*-(■,”)• (8) 


long bonds are generated in the transient states through 
which the system evolves. Fig. [G] shows the bond-length 
distribution for bonds of “shape” ( x , 0) in the simulations 
with the initial state of type Fig. |T])b) , at the point (in 
Fig. U) where the probability P(w x = 1) = 1/2. The 
results are compared to those in the VB solid ground 
state without domain walls. 

Clearly, the transient state with the decaying domain 
walls has a dramatically larger fraction of long bonds, 
which suggests a way to understand the short life time: 
The projected state evolves toward the ground state 
through transitional states that have a very high fraction 
of long bonds relative to the ground state dominated by 
very short bonds. Note that the initial state we are using 
only has short bonds, and by looking at the metastable 
domain wall quasi-eigenstate we also find that it is dom¬ 
inated by short bonds, very similar to the ground state. 
Our results show the presence of excited states which are 
very different and which can be components of the initial 
state (which they have to be in order to appear as a re¬ 
sult of the PQMC procedure), at least in part because of 
the over-completeness and non-orthogonality of the VB 
basis. The total fraction of bonds longer than Lj 4, i.e., 
those that can generate changes in the winding number, 
is graphed versus the system size in Fig. [7] It decays as 
L~ 2 in the case considered, while the decay is exponential 
in the ground state. 


IV. EFFECTIVE TWO-STATE MODEL AND 
REAL-TIME EVOLUTION 

It is not easy to relate our findings above to real-time 
evolution, where an isolated system would not decay to 
its ground state but go through a thermalization process 
with conserved excess energy density oc 1 /L. Also in this 
case one can define a time scale on which states with 
winding numbers different from the initial one become 


Starting with the initial state |V’( r )) = | t) we compute 
the probability P(t) = |(t \’ t P(, T ))\ 2 of the system staying 
in this initial state after time evolution using U(r) in 
Eq. ©. Defining the life time by P = 1/2 as in the 
preceding section gives the exact result 


2rv^+y2 rV 1 + r 2 - r 2 

vr+^+i ’ 


(9) 


where r = x/e. To leading order in r this becomes 


e 


—2re _ 


X 

2e 


( 10 ) 


Using the scaling behaviors found using the PQMC sim¬ 
ulations, e ~ L and T \/2 ~ L a , we must then have 


x 


Le~ L 


l + a 


(11) 


Going to real time, we instead (but equivalently) solve 
for P(i), obtaining 

P (I) = x2 X +e 2 sin (We 2 +£ 2 ) ■ (12) 

Here the oscillatory behavior is clearly different from 
what one would obtain in an infinite many-body system, 
where the initial state decays into many other states and 
no periodicity is expected for a thermalizing system in 
the thermodynamic limit. Nevertheless, one can define 
a rate of depletion of the initial state as the maximum 
P(i) divided by the time taken to reach this maximum. 
This gives the rate (for small x/e) v = 2 x 2 /(i re), which 
with Eq. m and e ~ L becomes 

v ~ Le~ 2L + . (13) 


While the approximation of the evolution of the quantum 
many-body state by just a two-level system, where the 
two levels represent entire sectors of states (blocks in the 
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Hamiltonian matrix), cannot be expected to be quanti¬ 
tatively accurate, the above simple calculation neverthe¬ 
less illustrates how the imaginary-time behavior can be 
qualitatively different from the corresponding real-time 
dynamics. The exponential decay rate obtained above is 
most likely correct, though details such as the power a 
in Eq. (1131) may not necessarily be accurate. It would be 
interesting to investigate this issue further by considering 
more sophisticated models with more than two states. 


V. RESULTS AT THE CRITICAL POINT 

Our studies of the life time of the winding number in 
the preceding sections, and the associated convergence 
of excited energies with system size to values consistent 
with the domain-wall energy per unit length, show un¬ 
ambiguously that the winding number is an emergent 
conserved quantity in the VB solid state. This in itself is 
perhaps not surprising, but our calculations have demon¬ 
strated the nature of the mechanism causing the winding 
transitions (topological fluctuations) and quantified the 
life time. An important question now is how the mecha¬ 
nism of winding number decay evolves as one approaches 
a critical point at which the VB solid order vanishes. 
Such a critical point can be reached in the J-Q model. 
In the present variant with six-spin Q-interaction (j5j) the 
critical value of the ratio q = Q/{ J + Q) is q c ~ 0.600^ 

In the theory of deconfined quantum-criticality^ there 
are two relevant diverging length scales upon approaching 
the critical point at q = q c . In addition to the standard 
correlation length £ ~ (q — q c )~ u , which can be defined, 
e.g., using the distance dependence of spin-spin or VB- 
VB correlations, there is a larger length scale £vb char¬ 
acterizing the thickness of domain walls inside the VB 
solid, with £vb ~ (q — q c )~ v and v' > The presence 
of two intrinsic physical length scales in the system makes 
the finite-size scaling of k very interesting, with dramatic 
cross-overs predicted^ Studying k in detail upon ap¬ 
proaching the critical point within the J-Q model is an 
interesting problem, which however is beyond the scope 
of the present study. We here just consider the wind¬ 
ing number conservation and associated critical domain 
wall energy as a function of the system size at the critical 
point. Critical domain walls have previously been studied 
in the classical dimer model as well as in spin-liquid wave 
functions defined using short VBs i 50 i 51 In those cases the 
winding number is conserved by construction. 

Figure [5] shows the excitation energy (critical domain 
wall energy) per unit length at q = 0.6, which is within 
error bars of the best known value of q c for this modeL 14 
PQMC results for w x = 1, computed in the manner dis¬ 
cussed in Sec. uni are shown versus 1 /L for three different 
cases of /3(L) = aL. For reference, results deep inside the 
VB solid phase (q = 1) are also shown here. They con¬ 
verge to a non-zero constant as L —> oo, with two choices 
of (3(L) seen to produce the same result. At q c , going 
from f3 = L/8 to L/ 4, we can see a clear decrease in 



FIG. 8 . (Color online) VB domain-wall energy per unit length 
as a function of the inverse system size graphed on a log-log 
scale. In the case of the strongly-ordered VB solid (5 = 1 ), 
the energy computed with two different inverse temperatures 
P(L) converges to the same non-zero value as L —> 00 , while 
at the critical point (q = 0.6 ss q c ) convergence of the energy 
for all L with increasing /3 is demonstrated (the same also 
holds true at q = 1 if still larger j3 is used). The converged 
energy decays as a power-law form ~ L~ b . The fitted line 
shown here has slope b = 1.80 ± 0 . 01 . 


the energy, while upon further reducing the temperature 
to f3 = 3L/8 there are no significant changes for any L 
within the error bars (which now are large for large sys¬ 
tems, making calculations at still higher (3 prohibitively 
expensive). The results shown were obtained with an ini¬ 
tial state of the type in Fig. [Q)b), but converged results 
from the type-(c) state are the same within error bars. 

Interestingly, a very good power-law behavior is seen 
at q C: k ~ L~ b with b ss 1.8, which corresponds to the 
(quasi) eigen-energy E 0 (w x = 1) ~ L x ~ b « L~ 0 & . The 
lowest singlet energy in the J-Q model at q c scales as 
1/L^ as expected for a critical point with dynamic ex¬ 
ponent z = 1. Thus, the critical domain wall energy is 
only slightly above the lowest singlet (and it should be 
noted here again that the momentum of all the states we 
are computing here is 0 , as in the ground state). 

The domain wall energy per unit length of a VB solid 
can be expressed as 


k = 


K 

X’ 


(14) 


where I\ is a stiffness constant describing the energy cost 
of a twist of the VB order parameter and A is the width 
of the region over which this twist is distributed. In the 
theory of deconfined quantum-criticality^ the VB stiff¬ 
ness in the thermodynamic limit scales as K ~ £ _1 ~ 
(q ~ Qc)-" upon approaching the critical point, while 
A must saturate at the domain-wall thickness discussed 
above. Thus, in systems with domain walls imposed 
through winding numbers, one can expect that 


1 1 


(.q-qcY+r 


(15) 



2 


scales as L) becomes 



FIG. 9. (Color online) Size dependence of the imaginary-time 
life time at q c of states with winding numbers w x = 1 and 
w x = L/2, using initial states of the type in Fig. [2b) and its 
Wx > 1 generalization. The lines drawn through the large-/, 
points have slopes a = 1.2 and a = 0.75 for w x = 1 and 
w x = //2, respectively. 


v ~ Le~ 2Ll ~ b+a . (17) 

At w x = 1 the exponent 1 — b ± a remains positive with 
the values of a and b obtained above, and, thus, the tran¬ 
sition rate out of the original w x sector in real time is ex¬ 
ponentially small. Interesting, for w x = L/2 our results 
suggest 1 — b + a ~ 0 (within estimated error bars) and if 
this exponent indeed vanishes, and if the two-state model 
is to be taken seriously (which perhaps is asking too much 
of it in this extreme case), then the high-winding states 
would actually have a shorter than exponentially long life 
time at the critical point. 

Going to coupling ratios q < q c , entering the antifer¬ 
romagnetic state of the J-Q model, the winding number 
sectors should become completely unstable. Indeed, here 
the energy above the ground state computed with differ¬ 
ent /3 = aL quickly decays to 0 and, unlike the results at 
q c in Fig. [51 it is not possible to discern any converged 
functional form. 


In standard finite-size scaling procedures at a critical 
point, 53 to relate the behavior of a quantity in the ther¬ 
modynamic limit as a critical point is approached to the 
behavior as a function of the system size exactly at the 
critical point, one simply replaces the correlation length 
by the system length L. In the present case, we can ar¬ 
gue that it is £vb that should be replaced by L , since 
this length-scale is the one reaching L first when q c is 
approached for finite L. We then obtain 

k ( 9 c ) ~ L~^l v '\ (16) 

and, therefore, with the exponent b defined in the anal¬ 
ysis of our results above (shown in Fig. [5]) we have 
6=1 + vjv'. Thus, we have extracted a rather pre¬ 
cise estimate of the exponent ratio v/v' ~ 0.80 ±0.01, 
where the error bar is one standard deviation of the slope 
of the fitted line in Fig. [5] (and we estimate that the error 
due to very small deviations of q = 0.6 from the true q c 
is smaller than the quoted statistical error). The only 
other estimate of this exponent ratio that we are aware 
of is v'/v = 1.20 ± 0.05, or v/v' = 0.83 ± 0.04, from an 
analysis of the emergent 1/(1) symmetry of the VB or¬ 
der parameter^ It is gratifying that these two estimates 
obtained in completely different ways are fully consistent 
with each other. 

Now consider the life time in imaginary time, calcu¬ 
lated as in Sec. IIIII Fig.[9]shows results for w x = 1 as well 
as the extreme case of w x = L/2. We find ~ L 1,2 for 
w x = 1 and Ti /2 ~ L 0 ' 75 for w x = L/2 (with error bars 
on the exponents of about 10%). In the effective real¬ 
time model, with the critical domain wall energy scaling 
as Ew — Ew -1 ~ kL ~ L 1-b , the form (fl3l) of the de¬ 
cay rate derived in the VBS (where the energy difference 


VI. SUMMARY AND DISCUSSION 

We have demonstrated a mechanism of decay of do¬ 
main walls in a VB solid state of a quantum spin system, 
through fluctuations of the topological winding number 
which effectively counts the number of domain walls in 
a system. The mechanism requires the wave function to 
contain VBs of length proportional to the system size, 
which is explicitly excluded in an effective description of 
VB solids with a QDM. The domain walls become sta¬ 
ble in the thermodynamic limit, or, in other words, the 
winding number is an emergent conserved quantum num¬ 
ber. The life time in imaginary time scales as a power of 
the system size and we have argued, based on a simple 
two-state model for two winding sectors, that this trans¬ 
lates into an exponentially small transition rate out of an 
initial winding-number sector in real time. 

At a critical point separating the VB solid and an 
antiferromagnetic ground state (the putative deconfmed 
quantum-critical point^), we have also found stable 
winding numbers in the thermodynamic limit and a 
power-law decay of the excitation energy with the sys¬ 
tem size. The energy decay exponent contains informa¬ 
tion on the spectrum of the effective 1 /( 1 ) gauge-field 
model describing this phase transition, and our quanti¬ 
tative results from the power-law scaling, along with ex¬ 
ponents obtained previously for other quantities^ lend 
support to the CP 1 field theory of deconfmed quantum- 
criticality^ It would be interesting to study in detail the 
domain wall energy in the whole range of states between 
the maximally ordered VB solid and the critical point, to 
make further quantitative comparisons with the theory. 

It would be interesting to investigate the consequences 
of long VBs in spin liquid phases as well (noting that the 
critical point we have consider corresponds to an alge¬ 
braic spin liquid at an isolated point). To our knowledge, 
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the expected quasi-degenerate topological multiplet has 
never been observed in SU(2) invariant Z 2 spin-liquid 
candidate Hamiltonians^— (which are amenable to VB 
descriptions and whose topological numbers should be 
similar to the even-odd winding discussed in Ref. [l9h . 
This could be due to the splitting being relatively large, 
perhaps still exponentially small but with a large pref¬ 
actor, or there could be some cross-over from a different 
form when the system size is still relatively small and the 
effects of longer VBs could be significant. It should be 
noted here that, because of the over-completeness of the 
VB basis, to fix the winding number of an initial state 
it must overlap with states also outside the subspace of 
the topological multiplet of a gaped spin liquid. Some of 
these states may also have an enhanced fraction of long 
VBs, as we found here in the transitional states out of 
simple short-bond domain wall states. Spin liquids in 
frustrated spin models cannot be investigated with the 


PQMC methods we have used here, due to sign prob¬ 
lems, but it may be possible to study this issue with exact 
diagonalization techniques based on VBs,— though the 
limitation in system size may make it difficult to clearly 
observe the effects of long bonds. 
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